#################################################################################### 
####################################################################################
####################################################################################
# Appendix G: Tables G.15-G.17, G.19, M.48, M.50-M.52
####################################################################################
####################################################################################
library(here)
library(kableExtra)
library(lmtest)
library(multiwayvcov)
####################################################################################
####################################################################################
## Read in survey data (from: 8_2020_analysis.R)
dat = readRDS(here("data/dat_unrestricted.rds"))
datb = readRDS(here("data/dat_restricted.rds"))

# Implement SWE of ATE: Demean all covariates, fully interact
cols = c("sex","age","education_level","Arua","Bushenyi","Ibanda","Jinja","Mbale","Mpigi","Pallisa","Mbarara","Shema",
         "village_distance_HF","village_distance_HOSP","village_distance_road","village_distance_electricity", "share2006")
dat[, paste0("c_", cols)] = lapply(dat[, cols], function(x) (x - mean(x))  )
covariates = c("treatment*c_village_distance_HF + treatment*c_village_distance_HOSP + treatment*c_village_distance_road + treatment*c_village_distance_electricity + 
                treatment*c_sex + c_age*treatment + c_education_level*treatment + treatment*c_share2006 + 
                c_Arua*treatment + c_Bushenyi*treatment + c_Ibanda*treatment + c_Jinja*treatment + c_Mbale*treatment + c_Mpigi*treatment + c_Pallisa*treatment + c_Mbarara*treatment",
               "c_Arua*treatment + c_Bushenyi*treatment + c_Ibanda*treatment + c_Jinja*treatment + c_Mbale*treatment + c_Mpigi*treatment + c_Pallisa*treatment + c_Mbarara*treatment")

# Implement SWE of ATE: Demean all covariates, fully interact
cols = c("sex","age","education_level","Arua","Ibanda","Pallisa","Shema",
         "village_distance_HF","village_distance_HOSP","village_distance_road","village_distance_electricity", "share2006")
datb[, paste0("c_", cols)] = lapply(datb[, cols], function(x) (x - mean(x))  )
covariatesb = c("treatment*c_village_distance_HF + treatment*c_village_distance_HOSP + treatment*c_village_distance_road + treatment*c_village_distance_electricity + 
                treatment*c_sex + c_age*treatment + c_education_level*treatment + treatment*c_share2006 + 
                c_Arua*treatment + c_Ibanda*treatment + c_Pallisa*treatment", 
                "c_Arua*treatment + c_Ibanda*treatment + c_Pallisa*treatment")

####################################################################################
## H1 & H2: Full Contact Work Index

# Unrestricted
for (i in names(dat[,c("engage_health_gov_index","engage_other_gov_index",
                       "engage_health_ngo_index","engage_other_ngo_index")])){ # Dependent Variables
  for (j in covariates) {
    model = paste(i,"~","treatment","+", j)
    
    # Run each model
    assign(x = paste("m",i,substr(j,1,1), sep = "."), 
           value = lm(as.formula(model), data = dat))
    # Output clustered SEs (county)
    assign(x = paste("c",i,substr(j,1,1),sep = "."), 
           value = coeftest(lm(as.formula(model), data = dat),
                            cluster.vcov(lm(as.formula(model), data = dat), dat$village_id)))
  }
}

# Restricted
for (i in names(datb[,c("engage_health_gov_index","engage_other_gov_index",
                        "engage_health_ngo_index","engage_other_ngo_index")])){ # Dependent Variables
  for (j in covariatesb) {
    model = paste(i,"~","treatment","+", j)
    
    # Run each model
    assign(x = paste("bm",i,substr(j,1,1), sep = "."), 
           value = lm(as.formula(model), data = datb))
    # Output clustered SEs (county)
    assign(x = paste("bc",i,substr(j,1,1),sep = "."), 
           value = coeftest(lm(as.formula(model), data = datb),
                            cluster.vcov(lm(as.formula(model), data = datb), datb$village_id)))
  }
}

####################################################################################
####################################################################################
## Table G.15

stargazer(m.engage_health_gov_index.t, bm.engage_health_gov_index.t, m.engage_other_gov_index.t, bm.engage_other_gov_index.t,  
          m.engage_health_ngo_index.t, bm.engage_health_ngo_index.t, m.engage_other_ngo_index.t, bm.engage_other_ngo_index.t, 
          
          se = list(c.engage_health_gov_index.t[,2], bc.engage_health_gov_index.t[,2], c.engage_other_gov_index.t[,2], bc.engage_other_gov_index.t[,2], 
                    c.engage_health_ngo_index.t[,2], bc.engage_health_ngo_index.t[,2], c.engage_other_ngo_index.t[,2], bc.engage_other_ngo_index.t[,2] ),
          
          p = list(c.engage_health_gov_index.t[,4], bc.engage_health_gov_index.t[,4], c.engage_other_gov_index.t[,4], bc.engage_other_gov_index.t[,4], 
                   c.engage_health_ngo_index.t[,4], bc.engage_health_ngo_index.t[,4], c.engage_other_ngo_index.t[,4], bc.engage_other_ngo_index.t[,4] ),
          
          keep = c("treatment$"), 
          
          order = c("$treatment$"), covariate.labels=c("Treatment"),
          
          type = "latex", out = "/home/jeremy/Dropbox/Dissertation/Manuscript/DCRT/tables/subset_controls2/engagement_1_c.tex",
          label = "tab:engagement_1_c", column.sep.width = "1pt", table.placement = "!ht",
          keep.stat = c("n"), dep.var.labels.include = F, no.space = T, model.numbers = T,
          title = "Effect of LG CHP on Political Engagement by Policy Area",
          notes = "Standard errors are clustered at the village level. $^{*}$p$<$0.1; $^{**}$p$<$0.05; $^{***}$p$<$0.01", dep.var.caption = "", table.placement = "!ht",notes.align = "l", notes.append = F, notes.label = "",
          column.labels = c("Govt Health","Govt General","NGO Health","NGO General"), column.separate = c(2,2,2,2),
          add.lines = list(c("Restricted", "No", "Yes", "No", "Yes", "No", "Yes", "No", "Yes")))

####################################################################################
####################################################################################
## Table M.48

stargazer(m.engage_health_gov_index.c, bm.engage_health_gov_index.c, m.engage_other_gov_index.c, bm.engage_other_gov_index.c,  
          m.engage_health_ngo_index.c, bm.engage_health_ngo_index.c, m.engage_other_ngo_index.c, bm.engage_other_ngo_index.c, 
          
          se = list(c.engage_health_gov_index.c[,2], bc.engage_health_gov_index.c[,2], c.engage_other_gov_index.c[,2], bc.engage_other_gov_index.c[,2], 
                    c.engage_health_ngo_index.c[,2], bc.engage_health_ngo_index.c[,2], c.engage_other_ngo_index.c[,2], bc.engage_other_ngo_index.c[,2] ),
          
          p = list(c.engage_health_gov_index.c[,4], bc.engage_health_gov_index.c[,4], c.engage_other_gov_index.c[,4], bc.engage_other_gov_index.c[,4], 
                   c.engage_health_ngo_index.c[,4], bc.engage_health_ngo_index.c[,4], c.engage_other_ngo_index.c[,4], bc.engage_other_ngo_index.c[,4] ),
          
          keep = c("treatment$"), 
          
          order = c("$treatment$"), covariate.labels=c("Treatment"),
          
          type = "latex", out = "/home/jeremy/Dropbox/Dissertation/Manuscript/DCRT/tables/subset_controls2/engagement_1_nc.tex",
          label = "tab:engagement_1_nc", column.sep.width = "1pt", table.placement = "!ht",
          keep.stat = c("n"), dep.var.labels.include = F, no.space = T, model.numbers = T,
          title = "Effect of LG CHP on Political Engagement by Policy Area (No covariates)",
          notes = "Standard errors are clustered at the village level. $^{*}$p$<$0.1; $^{**}$p$<$0.05; $^{***}$p$<$0.01", dep.var.caption = "", table.placement = "!ht",notes.align = "l", notes.append = F, notes.label = "",
          column.labels = c("Govt Health","Govt General","NGO Health","NGO General"), column.separate = c(2,2,2,2),
          add.lines = list(c("Restricted", "No", "Yes", "No", "Yes", "No", "Yes", "No", "Yes")))



rm(list = grep("^(c|m)\\.", ls(), value = T))
rm(list = grep("^(bc|bm)\\.", ls(), value = T))

####################################################################################
## Contact Work
## H1 & H2: Full Contact Work Index 

# Unrestricted
for (i in names(dat[,c("contact_health_gov_index","contact_health_mp_st", "contact_health_dc_st", "contact_health_lc_st", 
                       "contact_health_dho_st", "contact_health_gak_st", "contact_health_ngo_st")])){ # Dependent Variables
  for (j in covariates) {
    model = paste(i,"~","treatment","+", j)
    
    # Run each model
    assign(x = paste("m",i,substr(j,1,1), sep = "."), 
           value = lm(as.formula(model), data = dat))
    # Output clustered SEs (county)
    assign(x = paste("c",i,substr(j,1,1),sep = "."), 
           value = coeftest(lm(as.formula(model), data = dat),
                            cluster.vcov(lm(as.formula(model), data = dat), dat$village_id)))
  }
}

# Restricted
for (i in names(datb[,c("contact_health_gov_index","contact_health_mp_st", "contact_health_dc_st", "contact_health_lc_st", 
                        "contact_health_dho_st", "contact_health_gak_st", "contact_health_ngo_st")])){ # Dependent Variables
  for (j in covariatesb) {
    model = paste(i,"~","treatment","+", j)
    
    # Run each model
    assign(x = paste("bm",i,substr(j,1,1), sep = "."), 
           value = lm(as.formula(model), data = datb))
    # Output clustered SEs (county)
    assign(x = paste("bc",i,substr(j,1,1),sep = "."), 
           value = coeftest(lm(as.formula(model), data = datb),
                            cluster.vcov(lm(as.formula(model), data = datb), datb$village_id)))
  }
}

####################################################################################
####################################################################################
## Table G.16

stargazer(m.contact_health_gov_index.t, m.contact_health_lc_st.t, m.contact_health_dc_st.t, m.contact_health_mp_st.t, 
          m.contact_health_dho_st.t,m.contact_health_gak_st.t, m.contact_health_ngo_st.t, 
          
          se = list(c.contact_health_gov_index.t[,2], c.contact_health_lc_st.t[,2], c.contact_health_dc_st.t[,2], 
                    c.contact_health_mp_st.t[,2], c.contact_health_dho_st.t[,2],  
                    c.contact_health_gak_st.t[,2], c.contact_health_ngo_st.t[,2] ),
          
          p = list(c.contact_health_gov_index.t[,4], c.contact_health_lc_st.t[,4], c.contact_health_dc_st.t[,4], 
                   c.contact_health_mp_st.t[,4], c.contact_health_dho_st.t[,4],  
                   c.contact_health_gak_st.t[,4], c.contact_health_ngo_st.t[,4] ),
          
          keep = c("treatment$"), 
          
          order = c("$treatment$"), covariate.labels=c("Treatment"),
          
          type = "latex", out = "/home/jeremy/Dropbox/Dissertation/Manuscript/DCRT/tables/subset_controls2/engagement_2_urc.tex",
          label = "tab:engagement_2_urc", column.sep.width = "1pt", dep.var.caption = "", table.placement = "!ht",
          keep.stat = c("n"), dep.var.labels.include = F, no.space = T, model.numbers = T,
          title = "Effect of LG CHP on Health-related Contact by Political Actor",
          notes = "Standard errors are clustered at the village level. $^{*}$p$<$0.1; $^{**}$p$<$0.05; $^{***}$p$<$0.01", notes.align = "l", notes.append = F, notes.label = "",
          column.labels = c("Govt","Local", "District","MP", "District","National",
                            "NGO\\\\ & Index & Councilors & Chair & & Agency & Agency"))

####################################################################################
####################################################################################
## Table M.50

stargazer(m.contact_health_gov_index.c, m.contact_health_lc_st.c, m.contact_health_dc_st.c, m.contact_health_mp_st.c, 
          m.contact_health_dho_st.c,m.contact_health_gak_st.c, m.contact_health_ngo_st.c, 
          
          se = list(c.contact_health_gov_index.c[,2], c.contact_health_lc_st.c[,2], c.contact_health_dc_st.c[,2], 
                    c.contact_health_mp_st.c[,2], c.contact_health_dho_st.c[,2],  
                    c.contact_health_gak_st.c[,2], c.contact_health_ngo_st.c[,2] ),
          
          p = list(c.contact_health_gov_index.c[,4], c.contact_health_lc_st.c[,4], c.contact_health_dc_st.c[,4], 
                   c.contact_health_mp_st.c[,4], c.contact_health_dho_st.c[,4],  
                   c.contact_health_gak_st.c[,4], c.contact_health_ngo_st.c[,4] ),
          
          keep = c("treatment$"), 
          
          order = c("$treatment$"), covariate.labels=c("Treatment"),
          
          type = "latex", out = "/home/jeremy/Dropbox/Dissertation/Manuscript/DCRT/tables/subset_controls2/engagement_2_urnc.tex",
          label = "tab:engagement_2_urnc", column.sep.width = "1pt", dep.var.caption = "", table.placement = "!ht",
          keep.stat = c("n"), dep.var.labels.include = F, no.space = T, model.numbers = T,
          title = "Effect of LG CHP on Health-related Contact by Political Actor (No covariates)",
          notes = "Standard errors are clustered at the village level. $^{*}$p$<$0.1; $^{**}$p$<$0.05; $^{***}$p$<$0.01", notes.align = "l", notes.append = F, notes.label = "",
          column.labels = c("Govt","Local", "District","MP", "District","National",
                            "NGO\\\\ & Index & Councilors & Chair & & Agency & Agency"))

####################################################################################
####################################################################################
## Table G.17

stargazer(bm.contact_health_gov_index.t, bm.contact_health_lc_st.t, bm.contact_health_dc_st.t, bm.contact_health_mp_st.t, 
          bm.contact_health_dho_st.t,bm.contact_health_gak_st.t, bm.contact_health_ngo_st.t, 
          
          se = list(bc.contact_health_gov_index.t[,2], bc.contact_health_lc_st.t[,2], bc.contact_health_dc_st.t[,2], 
                    bc.contact_health_mp_st.t[,2], bc.contact_health_dho_st.t[,2],  
                    bc.contact_health_gak_st.t[,2], bc.contact_health_ngo_st.t[,2] ),
          
          p = list(bc.contact_health_gov_index.t[,4], bc.contact_health_lc_st.t[,4], bc.contact_health_dc_st.t[,4], 
                   bc.contact_health_mp_st.t[,4], bc.contact_health_dho_st.t[,4],  
                   bc.contact_health_gak_st.t[,4], bc.contact_health_ngo_st.t[,4] ),
          
          keep = c("treatment$"), 
          
          order = c("$treatment$"), covariate.labels=c("Treatment"),
          
          type = "latex", out = "/home/jeremy/Dropbox/Dissertation/Manuscript/DCRT/tables/subset_controls2/engagement_2_rc.tex",
          label = "tab:engagement_2_rc", column.sep.width = "1pt", dep.var.caption = "", table.placement = "!ht",
          keep.stat = c("n"), dep.var.labels.include = F, no.space = T, model.numbers = T,
          title = "Effect of LG CHP on Health-related Contact by Political Actor (Restricted)",
          notes = "Standard errors are clustered at the village level. $^{*}$p$<$0.1; $^{**}$p$<$0.05; $^{***}$p$<$0.01", notes.align = "l", notes.append = F, notes.label = "",
          column.labels = c("Govt","Local", "District","MP", "District","National",
                            "NGO\\\\ & Index & Councilors & Chair & & Agency & Agency"))


####################################################################################
####################################################################################
## Table M.51

stargazer(bm.contact_health_gov_index.c, bm.contact_health_lc_st.c, bm.contact_health_dc_st.c, bm.contact_health_mp_st.c, 
          bm.contact_health_dho_st.c,bm.contact_health_gak_st.c, bm.contact_health_ngo_st.c, 
          
          se = list(bc.contact_health_gov_index.c[,2], bc.contact_health_lc_st.c[,2], bc.contact_health_dc_st.c[,2], 
                    bc.contact_health_mp_st.c[,2], bc.contact_health_dho_st.c[,2],  
                    bc.contact_health_gak_st.c[,2], bc.contact_health_ngo_st.c[,2] ),
          
          p = list(bc.contact_health_gov_index.c[,4], bc.contact_health_lc_st.c[,4], bc.contact_health_dc_st.c[,4], 
                   bc.contact_health_mp_st.c[,4], bc.contact_health_dho_st.c[,4],  
                   bc.contact_health_gak_st.c[,4], bc.contact_health_ngo_st.c[,4] ),
          
          keep = c("treatment$"), 
          
          order = c("$treatment$"), covariate.labels=c("Treatment"),
          
          type = "latex", out = "/home/jeremy/Dropbox/Dissertation/Manuscript/DCRT/tables/subset_controls2/engagement_2_rnc.tex",
          label = "tab:engagement_2_rnc", column.sep.width = "1pt", dep.var.caption = "", table.placement = "!ht",
          keep.stat = c("n"), dep.var.labels.include = F, no.space = T, model.numbers = T,
          title = "Effect of LG CHP on Health-related Contact by Political Actor (Restricted; No covariates)",
          notes = "Standard errors are clustered at the village level. $^{*}$p$<$0.1; $^{**}$p$<$0.05; $^{***}$p$<$0.01", notes.align = "l", notes.append = F, notes.label = "",
          column.labels = c("Govt","Local", "District","MP", "District","National",
                            "NGO\\\\ & Index & Councilors & Chair & & Agency & Agency"))

rm(list = grep("^(c|m)\\.", ls(), value = T))
rm(list = grep("^(bc|bm)\\.", ls(), value = T))

####################################################################################
## H1 Secondary

# Unrestricted
for (i in names(dat[,c("contentious_index","socialcap_index",
                       "info_consumption_index","knowledge_score")])){ # Dependent Variables
  for (j in covariates) {
    model = paste(i,"~","treatment","+", j)
    
    # Run each model
    assign(x = paste("m",i,substr(j,1,1), sep = "."), 
           value = lm(as.formula(model), data = dat))
    # Output clustered SEs (county)
    assign(x = paste("c",i,substr(j,1,1),sep = "."), 
           value = coeftest(lm(as.formula(model), data = dat),
                            cluster.vcov(lm(as.formula(model), data = dat), dat$village_id)))
  }
}

# Restricted
for (i in names(datb[,c("contentious_index","socialcap_index",
                        "info_consumption_index","knowledge_score")])){ # Dependent Variables
  for (j in covariatesb) {
    model = paste(i,"~","treatment","+", j)
    
    # Run each model
    assign(x = paste("bm",i,substr(j,1,1), sep = "."), 
           value = lm(as.formula(model), data = datb))
    # Output clustered SEs (county)
    assign(x = paste("bc",i,substr(j,1,1),sep = "."), 
           value = coeftest(lm(as.formula(model), data = datb),
                            cluster.vcov(lm(as.formula(model), data = datb), datb$village_id)))
  }
}

####################################################################################
####################################################################################
## Table G.19

stargazer(m.contentious_index.t, bm.contentious_index.t, m.info_consumption_index.t, bm.info_consumption_index.t,  
          m.knowledge_score.t, bm.knowledge_score.t, m.socialcap_index.t, bm.socialcap_index.t, 
          
          se = list(c.contentious_index.t[,2], bc.contentious_index.t[,2], c.info_consumption_index.t[,2], bc.info_consumption_index.t[,2], 
                    c.knowledge_score.t[,2], bc.knowledge_score.t[,2], c.socialcap_index.t[,2], bc.socialcap_index.t[,2] ),
          
          p = list(c.contentious_index.t[,4], bc.contentious_index.t[,4], c.info_consumption_index.t[,4], bc.info_consumption_index.t[,4], 
                   c.knowledge_score.t[,4], bc.knowledge_score.t[,4], c.socialcap_index.t[,4], bc.socialcap_index.t[,4] ),
          
          keep = c("treatment$"), 
          
          order = c("$treatment$"), covariate.labels=c("Treatment"),
          
          type = "latex", out = "/home/jeremy/Dropbox/Dissertation/Manuscript/DCRT/tables/subset_controls2/engagement_3_c.tex",
          label = "tab:engagement_3_c", column.sep.width = "1pt", table.placement = "!ht",
          keep.stat = c("n"), dep.var.labels.include = F, no.space = T, model.numbers = T,
          title = "Effect of LG CHP on Secondary Measures of Political Engagement",
          notes = "Standard errors are clustered at the village level. $^{*}$p$<$0.1; $^{**}$p$<$0.05; $^{***}$p$<$0.01", dep.var.caption = "", table.placement = "!ht",notes.align = "l", notes.append = F, notes.label = "",
          column.labels = c("Contentious","Information","Knowledge","Membership"), column.separate = c(2,2,2,2),
          add.lines = list(c("Restricted", "No", "Yes", "No", "Yes", "No", "Yes", "No", "Yes")))

####################################################################################
####################################################################################
## Table M.52

stargazer(m.contentious_index.c, bm.contentious_index.c, m.info_consumption_index.c, bm.info_consumption_index.c,  
          m.knowledge_score.c, bm.knowledge_score.c, m.socialcap_index.c, bm.socialcap_index.c, 
          
          se = list(c.contentious_index.c[,2], bc.contentious_index.c[,2], c.info_consumption_index.c[,2], bc.info_consumption_index.c[,2], 
                    c.knowledge_score.c[,2], bc.knowledge_score.c[,2], c.socialcap_index.c[,2], bc.socialcap_index.c[,2] ),
          
          p = list(c.contentious_index.c[,4], bc.contentious_index.c[,4], c.info_consumption_index.c[,4], bc.info_consumption_index.c[,4], 
                   c.knowledge_score.c[,4], bc.knowledge_score.c[,4], c.socialcap_index.c[,4], bc.socialcap_index.c[,4] ),
          
          keep = c("treatment$"), 
          
          order = c("$treatment$"), covariate.labels=c("Treatment"),
          
          type = "latex", out = "/home/jeremy/Dropbox/Dissertation/Manuscript/DCRT/tables/subset_controls2/engagement_3_nc.tex",
          label = "tab:engagement_3_nc", column.sep.width = "1pt", table.placement = "!ht",
          keep.stat = c("n"), dep.var.labels.include = F, no.space = T, model.numbers = T,
          title = "Effect of LG CHP on Secondary Measures of Political Engagement (No covariates)",
          notes = "Standard errors are clustered at the village level. $^{*}$p$<$0.1; $^{**}$p$<$0.05; $^{***}$p$<$0.01", dep.var.caption = "", table.placement = "!ht",notes.align = "l", notes.append = F, notes.label = "",
          column.labels = c("Contentious","Information","Knowledge","Membership"), column.separate = c(2,2,2,2),
          add.lines = list(c("Restricted", "No", "Yes", "No", "Yes", "No", "Yes", "No", "Yes")))

rm(list = grep("^(c|m)\\.", ls(), value = T))
rm(list = grep("^(bc|bm)\\.", ls(), value = T))
